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ABSTRACT 

Explicit algebraic stress models that are valid for three-dimensional turbulent flows in 
non-inertial frames are systematically derived from a hierarchy of second-order closure mod- 
els. This represents a generalization of the model derived by Pope [J. Fluid Mech. 72 , 331 
(1975)] who based his analysis on the Launder, Reece and Rodi model restricted to two- 
dimensional turbulent flows in an inertial frame. The relationship between the new models 
and traditional algebraic stress models - as well as anistropic eddy viscosity models - is 
theoretically established. The need for regularization is demonstrated in an effort to explain 
why traditional algebraic stress models have failed in complex flows. It is also shown that 
these explicit algebraic stress models can shed new light on what second-order closure models 
predict for the equilibrium states of homogeneous turbulent flows and can serve as a useful 
alternative in practical computations. 


"This research was supported by the National Aeronautics and Space Administration under NASA Con- 
tract Nos. NAS 1-18605 and NAS 1-19480 while the second author was in residence at the Institute for 
Computer Applications in Science and Engineering (ICASE), NASA Langley Research Center, Hampton, 
VA 23681. 
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1. INTRODUCTION 


Turbulent flows contain length and time scales that can change dramatically from one 
flow configuration to the next. Consequently, two-equation turbulence models - wherein 
transport equations are solved for two turbulent fields that are directly related to the length 
and time scales - represent the simplest level of Reynolds stress closure that can be geometry 
independent. This explains the current popularity of the K — e model which, unlike the older 
mixing length models, does not require the specification of empirical length or time scales 
that must be adjusted in an ad hoc fashion from one flow to the next (see Launder and 
vSpalding 1974). However, despite this positive feature, the K — e model shares a common 
deficiency with the older mixing length models: it is based on the Boussinesq eddy viscosity 
hypothesis. It is well-known that eddy viscosity models are unable to properly describe 
turbulent flows with body force effects arising from curvature or a system rotation. Likewise, 
their inaccurate prediction of normal Reynolds stress differences make eddy viscosity models 
incapable of describing secondary flows in non-circular ducts (c.f. Lumley 1978, Launder 1990 
and Speziale 1991). This has led to the development of nonlinear stress-strain relationships 
in turbulence modeling that transcend the Boussinesq hypothesis - an area of research that 
has received varying degrees of attention during the past three decades. 

Early work on the development of turbulence models with a nonlinear stress-strain rela- 
tion tended to be empirical in nature, relying on analogies with non-Newtonian flows (see 
Rivlin 1957, Lumley 1970 and Saffman 1977). The similarities between the laminar flow of a 
non-Newtonian fluid and the mean turbulent flow of a Newtonian fluid have long been rec- 
ognized (see Hinze 1975 for an interesting discussion of this point). In these older empirical 
models, which typically were obtained by simple tensor invariance arguments, the Reynolds 
stresses were taken to be nonlinear polynomial functions of the mean velocity gradients. In 
recent years, nonlinear Reynolds stress models of this type have been obtained within the 
context of two-equation turbulence modeling by more formal expansion techniques incorpo- 
rating, for example, the Direct Interaction Approximation (DIA) and the Renormalization 
Group (RNG) (see Yoshizawa 1984, Speziale 1987, Rubinstein and Barton 1990, and Yakliot 
et al. 1992). These models, which are characterized by an explicit relationship between the 
Reynolds stress tensor and the mean velocity gradients (and possibly their time derivatives) 
have come to be referred to as a anisotropic eddy viscosity models.” 

During the 1970’s, second-order closure models became popular in which closure was 
achieved based on the Reynolds stress transport equation where history and nonlocal effects 
are accounted for. By means of an equilibrium hypothesis, in which the Reynolds stress 
convection and transport terms were neglected, algebraic stress models were obtained from 
second-order closures (see Rodi 1976). In these models, the Reynolds stresses were related 
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implicitly to the mean velocity gradients. This provided the first formal means, based on 
a higher-order closure, to justify the extension of the Boussinesq hypothesis to incorporate 
nonlinearities in the mean velocity gradients. However, this type of algebraic stress model 
is cumbersome to implement in complex flows since the stress-strain relation is not explicit; 
numerical stiffness problems can result from the need for successive matrix inversions at each 
iteration. 

In an interesting paper that has to a large extent gone unnoticed, Pope (1975) developed 
a nonlinear Reynolds stress model by invoking the same equilibrium hypothesis as Rodi 
(1976). However, Pope (1975) actually presented a methodology for obtaining an explicit 
relation for the Reynolds stress tensor from the implicit algebraic equations that Rodi (1976) 
obtained from the Launder, Reece and Rodi (1975) model. This methodology, which leads 
to explicit algebraic stress models, is based on the use of integrity bases from linear algebra. 
Due to the complexity of the algebra, Pope (1975) was only able to obtain a solution for 
two-dimensional turbulent flows which he accomplished using the Launder, Reece and Rodi 
model. 

The main purpose of the present paper is to extend the results of Pope (1975) to three- 
dimensional turbulent flows in non-inertial frames starting from a more general hierarchy of 
second-order closure models — a task that is now feasible computationally via symbolic manip- 
ulation. The relationship between these new explicit algebraic stress models and anisotropic 
eddy viscosity models, as well as the older nonlinear Reynolds stress models, will be estab- 
lished in a systematic fashion. There appears to be considerable confusion in the turbulence 
literature concerning the relationship between these various types of nonlinear Reynolds 
stress models that needs to be clarified. It will be shown that the new explicit algebraic 
stress models represent the equilibrium Reynolds stress anisotropies predicted by second- 
order closures in homogeneous turbulent flows. Consequently, beyond their potential use in 
practical turbulence calculations, the explicit algebraic stress models to be derived herein 
can be used to explore the predictive capabilities of a hierarchy of second-order closures in 
homogeneous turbulence. These issues will be discussed in more detail in the sections to 
follow. 

2. THEORETICAL BACKGROUND 

The incompressible turbulent flow of a viscous fluid will be considered where the velocity 
v and kinematic pressure P are decomposed into ensemble mean and fluctuating parts as 
follows: 

v = v -f u, P = P + p. (I) 
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The Reynolds stress tensor r,j = u,uj is a solution of the transport equation (c.f. Hinze 
1975): 

^ - r, k ^- + «>„ - sn - 2Si m (e mti r, t + e mK T,,) + ©J + rV'r.j (2) 

Ut uXk Ox^ 

which is valid in an arbitrary non-inertial reference frame that can undergo a rotation with 
angular velocity fi, relative to an inertial framing. In (2), D/Dt = d/dt + v • V is the 
convective derivative, ei jk is the permutation tensor, v is the kinematic viscosity, and 


$ij = P 


' du, duj ' 
dx, dxi 


(3) 


/ dui duj\ 
\dx k dx k j 


(4) 


vfi = - Q — ( UiUjUk + puiSjk + puj6i k ) (5) 

dx k 

are, respectively, the pressure-strain correlation, the dissipation rate tensor, and the turbu- 
lent transport term. Turbulence models based on the Reynolds stress transport equation (2) 
are referred to as “second-order closures” since (2) is obtained by taking a second moment 
of the fluctuating Navier-Stokes equation. 

Homogeneous turbulent flows in equilibrium, as well as regions of inhomogeneous turbu- 
lent flows where there is a production-equals-dissipation equilibrium, satisfy the constraints: 


where 


Dt 


- 0 


T>1 + = 0 


bij — 


- 

2K 


( 6 ) 

(7) 

( 8 ) 


is the anisotropy tensor ( K = ^TTiTTi is the turbulent kinetic energy). In physical terms, 
this is an equilibrium for which convective and transport effects can be neglected; it is the 
basic equilibrium hypothesis used in the derivation of algebraic stress models. Although 
it constitutes an idealization, it is comforting to know that this equilibrium hypothesis 
is achievable in interesting physical limits that include homogeneous shear flow and the 
logarithmic region of an equilibrium turbulent boundary layer. 

It follows from (6) that 

7 ^ T ij _ Tjj D 7t /q\ 

Dt ~ K Dt 1 ’ 
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( 10 ) 


and, hence, by making use of the contraction of (2) and (7), we obtain 


D Tjj 


(V-e) 


K 


where V = —Tijdvijdx 3 is the turbulence production and e = \ea is the scalar turbulent 
dissipation rate. The substitution of (7) and (10) into (2) yields the following equilibrium 
form of the Reynolds stress transport equation: 

CP £^)~~ — T t ,k~r{ “ T jk “ ~ + 2 Cl m ^C 7 nkj^"ik T (11) 

A ox k dx k 

The dissipation rate tensor can be split into isotropic and deviatoric parts as follows 

2 

£{j — ~£&ij D^ij* ( 1 ^) 

By making use of (8) and (12), we can rearrange (11) into the alternative form 


where 


(V — e)b lJ = —-K S tJ — I<(b ik Sjk + bjkSik 

2 — _ 

~b mn S mn Sij) A \bi k {iO j k -f 2c m k] fhn ) 

~\~bjk(u){k 4" 2e m fc ! fi m )] T — n,j 

- 1 (dvj dvA _ 1 (dv_ j dvj\ 

‘ 2 ydxj dx t ) ’ U,J 2 dx, y 

n ij = < h, J — D^xj ■ 


( 13 ) 


(14) 

(15) 


In all of the commonly used second-order closure models, II , 3 is modeled in the general form 
(see Lumley 1978, Reynolds 1987 and Speziale 1991): 


II,j = eAj(b) + KM tJ ki(b) 


dvk 

dxi 


(16) 


For non-inertial frames, dvk/dxi is replaced with dvk/dxi + e m /fcf2 m in (16). The substitu- 
tion of (16) (with its non-inertial correction) into (13) yields a closed system of algebraic 
equations for the determination of the Reynolds stress anisotropy in terms of the mean ve- 
locity gradients. This constitutes the general form of algebraic stress models; the algebraic 
stress model of Rodi (1976) is obtained when the Launder, Reece and Rodi model for n,j 
is introduced into (13) and non-inertial effects are neglected. It is clear that these algebraic 
stress models are implicit in nature since the Reynolds stress tensor appears on both sides 
of the equation. In the next section we will examine how explicit relations can be obtained 
from (13). 
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3. EXPLICIT ALGEBRAIC STRESS MODELS 


When models for I I 4jr that are tensorially linear in the anisotropy tensor b %3 are considered, 
it is possible to obtain an explicit expression for b{j in terms of the mean velocity gradients. 
This can be accomplished by techniques from linear algebra that make use of integrity bases 
(namely, polynomial representations for isotropic tensor functions). Pope (1975) was the first 
to show this by basing his analysis on the Launder, Reece and Rodi model which is the most 
commonly used model for TT t y that is a linear function of bij. We will extend the analysis of 
Pope (1975) to a more general hierarchy of linear models for IT^ where non-inertial effects are 
accounted for. Furthermore, our analysis will be conducted for three-dimensional turbulent 
flows. Pope (1975) restricted his analysis to two-dimensional turbulent flows because of the 
severe complexity of the algebra - an obstacle that we will overcome by using Mathematica ™ 
(Wolfram 1988). 

It can be shown that the most general form of (16) that is tensorially linear in the 
anisotropy tensor b{ 3 is given by (c.f. Reynolds 1987 and Speziale 1991): 


where 


n 


v 


-Ciefy + C 2 I<S tJ + C 3 I< 


^bikSjk T bjkSik 


2 , ~n 

^ ?7i7i 



-\-C 4 A (bikWjk + bjkW ik) 


IF {j — T Smji^lm 


(17) 


(18) 


is the absolute vorticity tensor (namely, the vorticity tensor lo 13 relative to an inertial frame). 
The coefficients C\ — C 4 can be functions of the invariants of b tJ and can depend on V fe (of 
course, (16) requires that, at most, C\ is a linear function of V/s). The model of Launder, 
Reece and Rodi (1975), the model of Gibson and Launder (1978) and the linearized model 
of Speziale, Sarkar and Gatski (1991) - where the quadratic part of the slow pressure-strain 
correlation is neglected - are special cases of (17). For these models, we have the following 
coefficients: 


Launder, Reece and Rodi Model 

C\ = 3.0, C 2 = 0.8, C 3 = 1.75, C 4 = 1.31 


(19) 


Gibson and Launder Model 

C, = 3.6, C 2 = 0.8, C 3 = 1.2, C 4 = 1.2 (20) 
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Spcziale , Sarkar and Gatski Model 


C, = 3.4 + i.SV/e, c, = 0.8 - 1.3 1I'J\ C, = 1.25, C, = 0.40 


( 21 ) 


where II b = bijbij. It should be noted that for sufficiently small anisotropies, more compli- 
cated models can be approximated by the form (17). 

The direct substitution of (17) into (13) yields the equation: 


K = ~9 T i ( C‘ 


3 ' ' } 

(■ 


- f — - 2 — 

Sfj T ( (7,3 2) \ bjkS jk T bjkSik -b mn S mn &ij 


T ( ('4 — 2) bik (^jk T T bjk ^u>ik T ^ 2 Cm ^ r ^ m j | 


where 


=(K-*) 

I< 


£ 

If we introduce the dimensionless, rescaled variables: 

1 
2 


Stj = tt^( 2 ~ 


1 


W- = - C 4 ) 




+ 

C 3 -2 


f— 1 

VG - 2 / 


, o 

-mj* J t, m 


6ii 


r 4 I “ij 


then (19) reduces to the simpler form 

6- = -s*. - (b- k s- k + 


( 22 ) 

(23) 

(24) 


(25) 

(26) 

(27) 

(28) 


In matrix form, (28) can be written as 

b* = -S* - (b*S* + S*b* - -{b*S*}I) + b’W’ - W*b* (29) 

3 

where {•} denotes the trace and I denotes the unit tensor. This is a linear algebraic equation 
for the determination of b* in terms of S* and W*; the solution to (29) is of the general 
form 

b* = f(S*,W*). (30) 

Form invariance under an orthogonal coordinate transformation (with rotation tensor Q) 
requires that 

Qf(S*,W*)Q r = f(QS*Q r ,QW*Q r ) (31) 
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i.e., that f be an isotropic tensor function of its arguments. It can be shown that the 
satisfaction of (31) requires that 

b* = Y G {X) T (A) (32) 

A 

where is the integrity basis for functions of a symmetric and antisymmetric tensor and 
G< A > are scalar functions of the irreducible invariants of S* and W*. For the case under 
consideration, the integrity basis is given by (see Spencer 1971 and Pope 1975): 


T (1 ) = S* 

T ( 2 ) = s*W* - w*s* 

X ( 3 ) _ S *2 _ 1{S* 2 }I 

X (4) = w* 2 - |{W* 2 }I 
T( 5 ) — w*s* 2 — s* 2 W* 


T (6) = w* 2 S* + S*W* 2 - §{S*W* 2 }I 
T (7) - W*S*W* 2 - W* 2 S*W* 

T< 8 ) = S*w*s* 2 - S* 2 W*S* 

T< 9 ) = W* 2 S* 2 + S* 2 W* 2 - §{S* 2 W‘ 2 }I 
T< 10 ) — w*s* 2 w * 2 — W* 2 S* 2 W* 


The irreducible invariants of S* and W* are: 


,1 = {S* 2 } (34) 

rn = {W* 2 } (35) 

„ = {S* 3 } (36) 

, 4 = {S*W 2 } (3T) 

<1, = {S' 2 W* 2 }. (38) 

Hence, we have 

G* A) = 772,773, V 4 , rjs), A = 1, 2, ... 10. (39) 


Pope (1975) showed that for two-dimensional flows, only and T^, as well as 771 

and 7/2, are independent. Consequently, the calculations become much simpler. Here, we 
will obtain the three-dimensional solution using the same methodology that Pope used. 

The direct substitution of (32) into (28) yields 


YjG { x) t ( a) 

A 


- Y £iaT ( a) - Y g (A) [t (A) s* + s*t ( a) - |{t (a) s*}i 

A A ^ 

_ X ( A )W* + W*T^ A ^j 


(40) 


where we have made use of the fact that = S* . Since is an integrity basis we know 
that 


t ( a ) s * _|_ s*T ( a) - -{T (a) S*}I = Y ^A 7 T (7) 


(41) 
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T^W* - W*T< A > = J2 -A- yT (n,) (42) 

7 

(i.e., any polynomial in S* and W* can be expanded in the integrity bases; see Spencer 1971 
and Pope 1975). The 10 x 10 matrices H Xl and <A 7 are obtained by successive applications 
of the Cayley-Hamilton Theorem (the elements of these matrices are given in Appendix A). 
By substituting (41) and (42) into (40), we obtain an equation from which G W 

can be 

determined: 

G< A > = S\\ - £G™tf 7 A + £GMJ 7A . (43) 

7 7 

This is a 10 x 10 linear system of equations for the determination of the G' A ' which can be 
written in the matrix form 

AG = B (44) 

where the components of A, G and B are given by 


A , 7 — <A 7 7A 7 T <A 7 (45) 

G\ = G (A) (46) 

Bx = A a (47) 

(see Appendix A). The solution to (44) is given by 

G (A) = A^. (48) 

In order to obtain a closed form expression for G^ x \ it is necessary to analytically invert the 
matrix A. The manual inversion of matrix A is not feasible in light of its highly complex 
structure; however, such a mathematical computation can now be done symbolically by 
Mathematical M (Wolfram 1988). The resulting inverse matrix A -1 is quite complex, but as 
noted in (48), only the first column of its elements are needed. The resulting solution for 
G (A) is: 

G (6) = 9/D 
G (7) = 9/D 

G (8) = 9/D (49) 


GW = -i (6 - 37?! - 21 7/2 - 2 7/ 3 + 307/4 ) /D 
GW = - (3 + 3t/, - 67/2 + 2t / 3 + 67 / 4 ) /D 
G (3) = (6 - 3t 7 , - 12t/2 - 2t / 3 - 6 7/4) /D 
GW = —3 (3t?i + 27/3 + 67 / 4 ) /D 
G (5) = 9/D 


G (9 > = 18/D 
G (10) = 0 


where the denominator D is given by 


15 


B - 3 - -7/1 + 7/j - —7/2 - 87/47/2 + 3?/ 2 - 7/3 + -7/17/3 


2 '" ' " 2 
-2t/ 2 7/3 + 21 7/4 + 247/5 + 27/! 7/4 - 6t/ 2 7/4 


(50) 
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Unfortunately, this three-dimensional form is significantly more complex than the two- 
dimensional form derived by Pope (1975). Nevertheless, the result presented in (49) can be 
shown to reduce to Pope’s form in the two-dimensional limit where the invariants 773 and 7/4 
are zero and 775 becomes equal to Using this reduction and the redundancies in the 

integrity basis (T (4) • • • T^ 10 ) are linear combinations of T'( I) ,T (2) and T f3 ); see Pope 1975), 
the resulting expression for G ^ in the two-dimensional limit becomes 


(7 (1) = — 

3 


3 

- 2 771 

- 6772 

G {2) = 

3 


3 

- 277, 

- 6772 

G (3) = — 

6 


3 - 

to 

1 

6772' 

The anisotropy tensor b* is then given by 




b * = ~ 3 -2,?- 6k [ S ’ + (S ' W ‘ “ W ‘ S ‘> " 2 ( S>! ” 5 (S ' 2}1 )] ’ (W ) 

which can be shown to be identical to the expression in Pope (1975) after using (25)-(27) 
along with the constants (19) of the Launder, Reece and Rodi model. It should be noted 
that Pope used an alternative form for the third basis in (52) - constructed from the two 
and three dimensional Kronecker delta - that is directly proportional to our T^ 3 k 

It is now clear that the model derived by Pope (1975) (which is a special case of (52)) 
actually constitutes the explicit solution of the algebraic stress model of Rodi for b in terms 
of S and uj. Whereas Pope’s analysis is only valid for two-dimensional flows in an inertial 
frame, we have now succeeded - by means of the representation theorem embodied in (32) 
- to extend his result to three-dimensional flows in non-inertial frames. 

Some comments are in order concerning how the three-dimensional solution obtained 
herein compares with the result derived recently by Taulbee (1992). The three-dimensional 
result of Taulbee only contains the bases T^\T^,T^,T( 6 \ where the coefficients 
G^ x \ for the most part, only depend on the invariant 7/2- This simplification results from the 
neglect of the second strain-dependent term in the brackets on the r.h.s. of (22) (namely, 
it is obtained by setting C 3 = 2). We feel that this is a questionable approximation. Our 
research on homogeneous shear flow indicates that the coefficient C 3 — 2 is nonzero and must 
be nearly half as large as the coefficient C 4 — 2 in order to obtain a good description of shear 
flows. Furthermore, by neglecting this term in (22), the coefficients of the solution become 
virtually independent of the strain rate S{j and the model collapses to the standard K — £ 
model in the absence of rotational strains. This stands in serious contradiction of the exact 
two-dimensional result (52) as well as our new three-dimensional result. 
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4. COMPARISON WITH OTHER NONLINEAR REYNOLDS STRESS 
MODELS 


The explicit algebraic stress model that we derived in Eqs. (32) and (49) formally con- 
stitutes the equilibrium anisotropy tensor predicted by a hierarchy of second-order closure 
models for homogeneous turbulence. These second-order closures are typically solved in con- 
junction with a modeled transport for the turbulent dissipation rate that is of the general 


D ^_ c £ v c £ \c — (-T — 
Dt - Ce1 7T ~ C ‘ 2 K + Cs dx t { e T '> dx 3 


where C e i,C e2 and C £ are either constants, or functions of the invariants of b, that vary 
from one model to the next. For homogeneous turbulence, the transport term on the r.h.s. 
of (53) vanishes and the turbulent kinetic energy is a solution of the transport equation 



e. 


(54) 


From (53)-(54), it follows that Vje achieves the following equilibrium value (see Speziale 
1991) 

Z = (55) 

« <7.i - 1 

that will be used to determine g in homogeneous turbulence. Hence, the explicit algebraic 
stress model derived herein can be utilized to systematically explore the equilibrium states 
predicted by a hierarchy of second-order closures in homogeneous turbulent flows. 

As discussed earlier, the explicit algebraic stress model given by Eqs. (32) and (49) 


represents an exact solution to the general implicit algebraic stress model (22) that includes 
the model of Rodi (1976) as a special case. Consequently, by a direct analysis of this explicit 
model, new insights can be gained concerning the problems that the traditional algebraic 
stress models have had in the calculation of complex turbulent flows. The denominator of 


the coefficients given in (49) contains a sum of positive and negative terms which has 
the potential to become zero rendering singular behavior. For equilibrium homogeneous 
turbulence this will not happen since j?i . . . Tfs are constrained by the transport equations 
for r,j and e to yield a well behaved solution (see Speziale and Mac Giolla Mhuiris 1989). 
However, when an algebraic stress model is applied to complex non-equilibrium turbulent 
flows, where the underlying assumptions invoked in its derivation are no longer valid, singular 
behavior can result. This becomes a distinct possibility in turbulent flows that have localized 


strain rates that are large. When the traditional algebraic stress models of the form (22) are 
applied to turbulent flows with large strain rates where D can approach zero, we would expect 
the model to be ill-behaved. This problem has been experienced in the practical calculation of 
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complex turbulent flows where the Rodi algebraic stress model has at times failed to converge 
when solved iteratively (Demuren and Rodi 1984 and Demuren, private communication). 
Hence, it is clear that there is the need to regularize these explicit algebraic stress models. 
Such an approach has been used recently in the kinetic theory with considerable success (see 
Rosenau 1989). 

The model derived herein can be regularized by a Pade approximation (this approach 
was recently used by Yakhot et al. 1992 in the derivation of their RNG based dissipation 
rate equation). The main idea behind this approach can be easily illustrated for the two- 
dimensional model given in (52). For this case, the coefficient on the r.h.s. of (52) can be 
rewritten as 

= (56) 

3 — 277 a — 6 ? 7 2 3 - 2t ? 2 + 6< 2 V ' 

where 

1 = (S’*S’*)‘ /2 . C = W'wr,) 1 ' 2 - (57) 

It is clear from (56) that, for sufficiently large strain rates 77 , singularities can occur; on the 
other hand, the rotational strains do not cause any problems. In an equilibrium homogeneous 
turbulent flow, 77 is typically less than one. For example, in homogeneous shear flow (with 
shear rate S), the parameter SK/e achieves an equilibrium value of approximately 5. For 
the Launder, Reece and Rodi model, this corresponds to 77 « 0.2 (see Speziale and Mac 
Giolla Mhuiris 1989). Hence, we want to replace (56) with a regular function that, for 77 
sufficiently less than one, is approximately the same. This can be accomplished by a Fade 
approximation. For example, to the first order we take 


I j ~ 


1 - 


1 + T} 2 


(58) 


Substituting (58) into (56) then yields the expression 

3 ^ 3(1 + 1} 2 ) 

3 - 2771—6772 ~ 3 + 7] 2 + 6( 2 r] 2 + 6C 2 
which is regular. Hence, we obtain a regularized model that is well-behaved for all strain 
rates and constitutes a good approximation to the original model within the equilibrium 
range that it formally applies. For practical applications, the regularized model 


(59) 


b* = 


3(1 +r? 2 ) 


3 + r) 2 + 6 ( 2 77 2 + 6( 2 


S* + (S’W - W'S‘) - 2 (s* ! - ~{S- 2 }l) 


(60) 


should be used instead of the form (52) first derived by Pope (1975) for the special case 
of the Launder, Reece and Rodi model. For the three-dimensional case, a more systematic 
Pade approximation to (50) must be made. This leads to rather cumbersome expressions 
that require further simplification and are currently under investigation. 


11 



We will now show the relationship between the explicit algebraic stress model derived 
herein and existing anisotropic eddy viscosity models. It is clear that (32), with the coeffi- 
cients (49), is of the general mathematical form 

( dv k \* 


b*; = A 


o 


ijkl 


, dxi 


(61) 


where ( dv k /dxi )* = S kl + W kl and A ljkl is a fourth rank tensor that is a function of S* and 
W*. This constitutes the general form of anisotropic eddy viscosity models; the standard 
eddy viscosity hypothesis of Boussinesq is recovered in the limit as Aij k i becomes an isotropic 
tensor. 

Our explicit algebraic stress model for three-dimensional turbulent flows is tensorially 
quartic in the mean velocity gradients with coefficients that are ratios of polynomials in the 

invariants of S and W. We will make use of the fact that 

2 


' %] 


KSij + 2a, Kb* 


*? 


(62) 


where = [(3 2 — /((^7 3 — 2). The normalized anisotropy tensor b can then be expanded 

in a Taylor series. To the first order in the mean velocity gradients, we have 

(63) 




- 9*. 


from the algebraic stress model (32) and (49). The direct substitution of (63) into (62), after 
making use of (25), yields 

m = \ks„ - 2c;!As„ ( 64 ) 

where C* = \g (| — C<i). This is the eddy viscosity form of the standard K — e model (see 
Launder and Spalding 1974) with a coefficient C* that, for the second-order closures consid- 
ered herein, is close to the traditional value of 0.09. When quadratic terms are maintained, 
we have 

Kj = -s-j - ) + 2 (s^. - (<») 

which, in an inertial frame, yields 

2 „ K 3 -= - - I< 3 ^ 


( 66 ) 


Tij = —KSij — 2C* — Sij — fa —^-{Sik^kj + Sjk^ki) + fa 2~ (SiicSkj „ S mn S mn S, 

where fa = \g 2 ( 2 - C 4 ) (§ - C 2 ) and fa = g 2 { 2 - C 3 ) (| - C 2 ). This is identical in form 
to the nonlinear K — £ model of Speziale (1987) when convective effects are neglected. The 
nonlinear K — e model is given by 


K A 


{ i3 


= ^K6, i -2C^S ii -‘iC D C H ^ 


A ' 3 /i 


Sij ^ Smm &ij 


(67) 


—4CdC 2 -^~ (s ik S kj l ^S mn iJ m n Ll ij 




r S mri 
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where Co and are constants and 


£ _ DJxjj _ o 'vij, _ dvj_- ~ 

13 Dt dx-k *” k} dxk * kl 


( 68 ) 


is the Oldroyd derivative of S. When DSij/Dt is neglected, it is a simple matter to show 
that (67) is of the same form as (66). 

Eq. (66) is also of the same form as the anisotropic eddy viscosity models of Yoshizawa 
(1984) and Rubinstein and Barton (1990) with one exception: it does not contain a term of 
the type 


ft 


^ 2 k ^ kj ^ ^ mn^ nm ®ij 


)■ 


(69) 


This kind of term was eliminated by Speziale (1987) because it yields an erroneous prediction 
for isotropic turbulence subjected to a solid body rotation. More precisely, a term of the 
form (69) predicts that an initially isotropic turbulence develops Reynolds stress anisotropies 
when subjected to a solid body rotation - a result that is in contradiction of physical and 
numerical experiments (see Wigeland and Nagib 1978 and Speziale, Mansour and Rogallo 
1987). It is encouraging that this consistency constraint, concerning the vanishing of (69) 
to the second order in the mean velocity gradients, has now been obtained directly by a 
systematic analysis of the Reynolds stress transport equation. 

It is interesting to note that the quadratic approximation (66) to the three-dimensional 
algebraic stress model is of the same tensorial form as the full nonlinear two-dimensional 
model (60). However, due to the strain-dependent coefficients in (60), more physics is ac- 
counted for. For example, Eq. (60) properly predicts that in a rapidly rotating frame - where 
fi, and hence tends to infinity - &*• and V go to zero. This is the restabilization effect 
of turbulence in a rapidly rotating frame that the usual anisotropic eddy viscosity models 
are unable to predict (see Speziale, Gatski and Mac Giolla Mhuiris 1990). Hence, we feel 
that the two-dimensional algebraic stress model (60) should be tried in the future as an 
alternative to the more commonly used quadratic anisotropic eddy viscosity models. Some 
applications of (60) will be considered in the next section. 


5. ILLUSTRATIVE EXAMPLES 

In order to demonstrate the efficacy of the regularized algebraic stress model (60) derived 
in this study, we will consider a few applications to non-trivial, two-dimensional turbulent 
fiows involving shear and rotation. We will base our calculations on the Speziale, Sarkar and 
Gatski (1991) second-order closure (hereafter, referred to as the SSG model). The regularized 
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algebraic stress model can be rewritten in the form 


‘v 


2 j,r 

= 3 h 4 


*3 


6(1 +r 1 2 )a l K 
3 + 7? 2 + 6£V + 6C 2 L 


1 


S’, + + SJ^) - 2 [s-.Sl, - ^S- U S’ U 8„ 


(70) 

where «i = ( C 2 — |) / (C3 — 2) and S*, W*, •q and £ are as defined in (25), (26), and 
(57). Some remarks are needed concerning how V/e and II are evaluated in the model 
coefficients (21) for the SSG model. V/e is calculated using (55) which is formally valid 
for equilibrium homogeneous turbulent flows; we take lit, ~ 0.11 which is the universal 
equilibrium value predicted by the SSG model for two-dimensional homogeneous turbulence 
(see Speziale, Sarkar and Gatski 1991). This yields the following choice of constants 


Ci = 6.80, C 2 = 0.36, C 3 = 1.25, C 4 = 0.40, g = 0.233 


(71) 


for the explicit algebraic stress model corresponding to the SSG second-order closure. 

The explicit algebraic stress model (70) is solved in conjunction with modeled transport 
equations for K and e. These equations are of the same general form as those used in the 
K — £ model and are given by: 

DI< 


Dt 


dvi 

d t 

f l/T 

dlO 


(72) 

dxj 

£+ fe( 


dxij 

1 

dv t 

„ £ 2 

d 

( Vt 


(73) 

dxj 

c ' 2 k + 

dxi 

u 

dxi) 

1, cr e 

= 1.3, C t 

-1 = 

1.44 

and C c2 — 1-83. 

Hence, the 


— = -C —r 

Dt ° £l K ! - 


explicit algebraic stress model derived herein formally constitutes a two-equation turbulence 
model; the K — £ model is recovered when (70) is linearized with respect to the mean velocity 
gradients. While this explicit algebraic stress model is a two-equation model, it incorporates 
much more physics than the K — e model since it is consistent with second-order closures in 
the limit of equilibrium homogeneous turbulent flows. 

The first example that we will consider is the case of homogeneous shear flow in a rotating 
frame. An initially decaying isotropic turbulence is, at time t = 0, subjected to a uniform 
shear rate S in a reference frame rotating steadily with angular velocity fi. The corresponding 
mean velocity gradient tensor is given by 


dz\ 

dxj 


0 5 1 
= I 0 0 
0 0 


(74) 


and the angular velocity of the reference frame is given by 12, = (0,0,12) (see Figure 1). In 
Table 1, the equilibrium values predicted by the new explicit algebraic stress model (ASM) 
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for homogeneous shear flow are compared with the experimental data of Tavoularis and 
Corrsin (1981) as well as with the predictions of the SSG model and the standard K — £ 
model. From these results it is clear that: (a) the new explicit ASM performs far better than 
the standard K — e model, and (b) the new explicit ASM yields results in close proximity 
of the full SSG second-order closure model. The small differences between the new ASM 
model and the full SSG model arises from two sources: the neglect of the quadratic return 
to isotropy term and the Pade approximation that was implemented to regularize the model. 
However, it is clear from these results that the new regularized ASM model constitutes an 
excellent approximation to the full SSG model for equilibrium flows. 

Now, we will compare the predictions of these three models for the time-evolution of the 
turbulent kinetic energy (K* = K/K 0 , t* = St) with the large-eddy simulations (LES) of 
Bardina, Ferziger and Reynolds (1983) for rotating homogeneous shear flow. In Figures 2(a)- 
2(c), the model predictions for three rotation rates (f l/S = 0, $l/S = 0.5 and f 1/S = —0.5) 
are compared with the LES which is for an initial condition of £o/ SKq = 0.296. From these 
results it is clear that the new algebraic stress model does an excellent job in capturing the 
trends of the LES. It yields results that are far superior to the standard K — £ model and 
are in close proximity to the SSG model. The main discernible difference between the new 
algebraic stress model and the SSG second-order closure is during the early transient where 
it responds more abruptly to the application of the shear since it does not account directly 
for relaxation effects. However, it is important to note how the new ASM is far superior to 
the standard K — e model in responding to changes in the rotation rate. This is illustrated in 
Figures 3(a)-3(c) where it can be seen that the new algebraic stress model is able to capture 
the effect of rotations on homogeneous shear flow. In contrast to these results, the standard 
K — £ erroneously yields the same results for all rotation rates — a deficiency tied to the 
Boussinesq eddy viscosity hypothesis. 

The last example that we will consider is the case of fully-developed turbulent channel 
flow subjected to a spanwise rotation with constant angular velocity f! (see Figure 4). We will 
consider the experimental test case of Johnston, Halleen and Lezius (1972) for a Reynolds 
number Re = 11,500 and a rotation number Ro = 0.21 (here, Re = UoH/u and Ro = 
HH/Uo where Uo is the bulk mean velocity). In Figure 5(a) the prediction for the mean 
velocity profile obtained from the SSG second-order closure by Speziale, So and Younis 
(1992) is compared with the experimental data of Johnston, Halleen and Lezius (1972). The 
computations were done using law of the wall boundary conditions. In Figure 5(b), the same 
comparisons are made for the new explicit algebraic stress model (60) derived herein. It is 
clear from these results that the new ASM yields an asymmetric mean velocity profile that 
is very similar to that obtained from the full second-order closure. Due to the use of wall 
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functions - as well as the neglect of roll instabilities (see Speziale, So and Younis 1992) - the 
specific quantitative comparisons are not that good. The important point, however, is that 
the new algebraic stress model yields an asymmetric mean velocity profile in close proximity 
to that obtained from a full second-order closure. On the other hand, simpler two-equation 
models like the standard K — e model erroneously predict that the mean velocity profile is 
unaffected by a system rotation and remains symmetric (see Speziale 1991). The improved 
predictions of the new explicit ASM in rotating shear flows arises from the incorporation of 
rotational strains through the terms W* and (. 

6. CONCLUSION 

Explicit algebraic stress models for two and three dimensional turbulent flows, in uon- 
inertial frames, have been obtained for a hierarchy of second-order closure models that are 
tensorially linear in the Reynolds stress anisotropy. These models were obtained using the 
standard local equilibrium hypothesis and, therefore, constitute the explicit solution to the 
traditional algebraic stress models generalized to include non-inertial effects and a range of 
pressure-strain models. They also formally represent the equilibrium states predicted by this 
hierarchy of second-order closures in homogeneous turbulent flows. A direct examination of 
these explicit models has shed new light on the limitations of the traditional algebraic stress 
models in applications to complex turbulent flows. For localized strain rates that are large, 
traditional algebraic stress models can become singular. These models need to be regularized 
- a task that can be achieved by means of a Pade approximation as demonstrated in this 
paper. 

The results of this study have presented the first definitive evidence as to why the tradi- 
tional algebraic stress models are ill-behaved, yielding numerical problems in many applica- 
tions. These models should eventually be abandoned in favor of regularized versions of the 
explicit algebraic stress models derived herein. The full three-dimensional form of the explicit 
models is rather complicated and work is underway to simplify them by means of a ratio- 
nal approximation procedure. However, in the mean time, the regularized two-dimensional 
form (60) can be used in practical applications. This simplified model is formally valid for 
two-dimensional equilibrium flows and collapses to the recently proposed anisotropic eddy 
viscosity models in the limit of small strain rates. However, while its general tensorial form 
is as simple as the commonly used anisotropic eddy viscosity models, it incorporates much 
more physics since the coefficients depend nonlinearly on both rotational and irrotational 
strains. A version of (60) based on the SSG second-order closure was applied to rotating 
shear flows herein with encouraging results. These developments, when combined with some 
recent improvements in the modeling of the turbulent dissipation rate, can lead to a new 
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generation of two-equation models that can serve as a useful companion to second-order 
closures in the calculation of complex turbulent flows. 
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APPENDIX A 


The elements of the matrices a and J 7 A are given in terms of rj, as 
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The elements of the matrix A are a linear combination of the elements of the matrices 
H-yX and J-y\, as well as the identity matrix. The matrix structure is then given by 

— 1 Au Aj 3 A14 0 Aie An Ais A19 0 

A21 — 1 0 0 A25 A26 A27 A28 0 A210 

A 3 i 0 —1 0 A 35 A 36 0 A 38 A39 A310 

0 0 0 —1 A45 A 46 A 47 0 A49 A410 

0 A52 A53 0 — 1 0 A57 0 A59 A510 

0 0 0 A 64 0 — 1 A67 A 68 A69 0 

0 0 0 0 0 A 76 —1 0 0 A 7 io 

0 0 0 0 A 85 0 0 —1 0 A 8 io 

0 0 0 0 A 95 Age 0 0 1 A 9 10 

0 0 0 0 0 0 A 107 0 A 109 —I 

with the nonzero elements determined from 

A 7 a = - Hx-y + Jx-r ( A4 ) 

The vector B is given by 
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Equilibrium 

Values 

Standard 
K — e Model 

New ASM 
Model 

SSG Model 

Experimental 

Data 

bn 

0 

0.204 

0.218 

0.20 

b 12 

-0.217 

-0.157 

-0.163 

-0.15 

&22 

0 

-0.149 

-0.146 

-0.14 

b 33 

0 

-0.055 

-0.072 

-0.06 

SK/e 

4.82 

6.02 

5.76 

6.0 


Table 1. Comparison of the model predictions with the experimental equilibrium values in 
homogeneous shear flow measured by Tavoularis and Corrsin (1981). 
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Figure 2. Time evolution of the turbulent kinetic energy in rotating homogeneous shear 
flow: Comparison of the model predictions with the large-eddy simulations of Bardina ct nl. 
(1983). (a) n/S = 0, (b) Sl/S = 0.5 and (c) Sl/S = -0.5. 
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Figure 2. Time evolution of the turbulent kinetic energy in rotating homogeneous shear 
flow: Comparison of the model predictions with the large-eddy simulations of Bardina et al. 
(1983), (a) Cl/S = 0, (b) Cl/S = 0.5 and (c) Cl/S = —0.5. 
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Figure 2. Time evolution of the turbulent kinetic energy in rotating homogeneous shear 
flow: Comparison of the model predictions with the large-eddy simulations of Bardina cl al. 

(1983). (a) 0/5 = 0, (b) 0/5 = 0.5 and (c) 0/5 = -0.5. 
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Figure 3. Time evolution of the turbulent kinetic energy in rotating homogeneous shear How: 
(a) large-eddy simulations of Bardina et al. (1983), (b) new explicit algebraic stress model, 
and (c) standard K — £ model. 
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Figure 3. Time evolution of the turbulent kinetic energy in rotating homogeneous shear (low: 
(a) large-eddy simulations of Bardina ct nl. (1983), (b) new explicit algebraic stress model, 
and (c) standard K — e model. 
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Figure 3. Time evolution of the turbulent kinetic energy in rotating homogeneous shear flow: 
(a) large-eddy simulations of Bardina el al (1983), (b) new explicit algebraic stress model, 
and (c) standard K — e model. 
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I'igure 5. Comparison of the model predictions (-) for the mean velocity profile with the ex- 
perimental data (□) of Johnston et al. (1972) for rotating channel flow (Re = 11,500, Ro = 
0.21): (a) SSG second-order closure (taken from Spezialc et al. 1992) and (b) new explicit 
algebraic stress model. 
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